cimport numpy
import numpy

DTYPE = numpy.float
ctypedef numpy.float_t DTYPE_t

cpdef inline IntVec cross(IntVec a, IntVec b):
    return a.cross(b) 

cdef class IntVec:
    """
    This class represents a Vector in 3D space.
    """

    # Declared in the .pxd file.
    #cdef public long x, y, z
    
    def __init__(self, x=0, y=0, z=0):
        """Constructor for a IntVec."""
        self.x = x
        self.y = y
        self.z = z

    def __reduce__(self):
        """
        Implemented to facilitate pickling of the Vec extension type.
        """
        d = {}
        d['x'] = self.x
        d['y'] = self.y
        d['z'] = self.z

        return (IntVec, (), d)

    def __setstate__(self, d):
        self.x = d['x']
        self.y = d['y']
        self.z = d['z']

    def __str__(self):
        return '(%f, %f, %f)'%(self.x, self.y, self.z)

    def __repr__(self):
        return 'IntVec(%g, %g, %g)'%(self.x, self.y, self.z)

    def __add__(IntVec self, IntVec p):
        return IntVec_new(self.x + p.x, self.y + p.y, self.z + p.z)

    def __sub__(IntVec self, IntVec p):
        return IntVec_new(self.x - p.x, self.y - p.y, self.z - p.z)

    def __mul__(self, m):
        cdef IntVec rself
        cdef int rm
        if(isinstance(self, IntVec)):
            rself = self
            rm = m
            return rself.mul(rm)
        rself = m
        rm = self
        return rself.mul(rm)

    cpdef IntVec mul(self, m):
        return IntVec_new(self.x*m, self.y*m, self.z*m)
    
    cpdef richcmp(self,IntVec p,int oper):
        if oper == 2: # ==
            if self.x == p.x and self.y == p.y and self.z == p.z:
                return True
            return False
        elif oper == 3: # !=
            if self.x == p.x and self.y == p.y and self.z == p.z:
                return False 
            return True
        else:
            raise TypeError('No ordering is possible for Vecs.')
    
    cpdef IntVec iadd(self,IntVec p):
        self.x += p.x
        self.y += p.y
        self.z += p.z
        return self
    
    cpdef IntVec isub(self,IntVec p):
        self.x -= p.x
        self.y -= p.y
        self.z -= p.z
        return self
    
    cpdef IntVec imul(self,m):
        self.x *= m
        self.y *= m
        self.z *= m
        return self
    
    cpdef IntVec idiv(self,m):
        self.x /= m
        self.y /= m
        self.z /= m
        return self
    
    #cpdef double getitem(self,int i):
    #    return (&self.x)[i]
    
    #cpdef setitem(self,int i,val):
    #    (&self.x)[i] = val
    
    cpdef add_inplace(self, IntVec p):
        self.x += p.x
        self.y += p.y
        self.z += p.z
        
    cpdef IntVec add(self, IntVec p):
        return IntVec_new(self.x + p.x, self.y + p.y, self.z + p.z)
    
    def __div__(IntVec self, m):
        return self.mul(1/m)

    def __abs__(self):
        return self.length()

    def __neg__(self):
        return IntVec_new(-self.x, -self.y, -self.z)

    def __richcmp__(IntVec self, IntVec p, int oper):
        if oper == 2: # ==
            if self.x == p.x and self.y == p.y and self.z == p.z:
                return True
            return False
        elif oper == 3: # !=
            if self.x == p.x and self.y == p.y and self.z == p.z:
                return False 
            return True
        else:
            raise TypeError('No ordering is possible for Vecs.')

    def __iadd__(self, IntVec p):
        self.x += p.x
        self.y += p.y
        self.z += p.z
        return self

    def __isub__(self, IntVec p):
        self.x -= p.x
        self.y -= p.y
        self.z -= p.z
        return self

    def __imul__(self, m):
        self.x *= m
        self.y *= m
        self.z *= m
        return self

    def __idiv__(self, m):
        self.x /= m
        self.y /= m
        self.z /= m
        return self
        
    #def __getitem__(self, int i):
    #    return (&self.x)[i]
    
    #def __setitem__(self, int i, val):
    #    (&self.x)[i] = val

    cpdef set(self, x, y, z):
        """Set the position from a given array.
        """
        self.x = x
        self.y = y
        self.z = z

    cpdef numpy.ndarray asarray(self):
        """Return a numpy array with the coordinates."""
        cdef numpy.ndarray[DTYPE_t, ndim=1] r = numpy.empty(3)
        r[0] = self.x
        r[1] = self.y
        r[2] = self.z
        return r

    cpdef inline double norm(self):
        """Return the square of the Euclidean distance to this Vec."""
        return (self.x*self.x + self.y*self.y + self.z*self.z)

    cpdef double length(self):
        """Return the Euclidean distance to this Vec."""
        return sqrt(self.norm())

    cpdef double dot(self, IntVec p):
        """Return the dot product of this Vec with another."""
        return self.x*p.x + self.y*p.y + self.z*p.z

    cpdef IntVec cross(self, IntVec p):
        """Return the cross product of this Vec with another, i.e.
        `self` cross `p`."""
        return IntVec_new(self.y*p.z - self.z*p.y, 
                       self.z*p.x - self.x*p.z,
                       self.x*p.y - self.y*p.x)

    cpdef double distance(self, IntVec p):
        """Return the distance between this Vec and p"""
        return sqrt((p.x-self.x)*(p.x-self.x)+
                    (p.y-self.y)*(p.y-self.y)+
                    (p.z-self.z)*(p.z-self.z))
                    
cpdef inline IntVec IntVec_new(x=0, y=0, z=0):
    cdef IntVec p = IntVec.__new__(IntVec)
    p.x = x
    p.y = y
    p.z = z
    return p

cdef inline IntVec IntVec_sub(IntVec pa,IntVec pb):
    return IntVec_new(pa.x-pb.x, pa.y-pb.y, pa.z-pb.z)

cdef inline IntVec IntVec_add(IntVec pa,IntVec pb):
    return IntVec_new(pa.x+pb.x, pa.y+pb.y, pa.z+pb.z)


cdef inline double IntVec_length(IntVec p):
    return sqrt(p.x*p.x + p.y*p.y + p.z*p.z)

cdef inline double IntVec_length2(IntVec p):
    return p.x*p.x + p.y*p.y + p.z*p.z

cdef inline double IntVec_distance(IntVec pa,IntVec pb):
    return sqrt((pa.x-pb.x)*(pa.x-pb.x) +
                (pa.y-pb.y)*(pa.y-pb.y) + 
                (pa.z-pb.z)*(pa.z-pb.z)
                )

cdef inline double IntVec_distance2(IntVec pa,IntVec pb):
    return ((pa.x-pb.x)*(pa.x-pb.x) + (pa.y-pb.y)*(pa.y-pb.y) + 
                    (pa.z-pb.z)*(pa.z-pb.z))

